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Abstract. 

We perform a systematic mimerical study, based on the time-dependent Gross- 
Pitaevskii equation, of jet formation in collapsing and exploding Bose-Einstein 
condensates as in the experiment by Donley et al. [2001 Nature 412 295]. In the actual 
experiment, via a Feshbach resonance, the scattering length of atomic interaction was 
suddenly changed from positive to negative on a pre-formed condensate. Consequently, 
the condensate collapsed and ejected atoms via explosion. On a disruption of collapse 
by suddenly changing the scattering length to zero a radial jet of atoms was formed 
in the experiment. We present a satisfactory account of jet formation under the 
experimental conditions as well as make predictions beyond experimental conditions 
which can be verified in future experiments. 



PACS numbers: 03.75.-b, 03.75.Nt 
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1. Introduction 

Recent successful observation of Bose-Einstein condensates (BEC) of trapped alkali 
atoms has initiated the intensive study of different novel phenomena. On the theoretical 
front, numerical simulation based on the time-dependent nonlinear mean-field Gross- 
Pitaevskii (GP) equation [1] is well under control and has provided a satisfactory account 
of some of these phenomena. 

Since the detection of BEC of ^Li atoms with attractive interaction, one problem of 
extreme interest is the dynamical study of the formation and decay of BEC for attractive 
atomic interaction [2]. In general a attractive condensate larger than a critical size is 
not dynamically stable [2]. However, if such a condensate is "prepared" or somehow 
made to exist it experiences a dramatic collapse and explodes emitting atoms. 

A dynamical study of the collapse has been performed by Donley et al. [3] on an 
attractive ^^Rb BEC [4] in an axially symmetric trap, where they manipulated the inter- 
atomic interaction by changing the external magnetic field exploiting a nearby Feshbach 
resonance [5] . In the vicinity of a Feshbach resonance the atomic scattering length a can 
be varied over a huge range by adjusting the external magnetic field. Consequently, they 
changed the sign of the scattering length, thus transforming a repulsive condensate of 
^^Rb atoms into an attractive one which naturally evolves into a collapsing and exploding 
condensate. Donley et al. have provided a quantitative estimate of the explosion of this 
BEC by measuring different properties of the exploding condensate. 

It has been realized that many features of the experiment by Donley et al. [3] 
can be described by the mean- field GP equation [6-17]. However, we are fully aware 
that there are features of this experiment which are expected to be beyond mean-field 
description. Among these are the distribution of number and energy of emitted high- 
energy (~ 10~^ Kelvin) uncondensed burst atoms reported in the experiment. Although 
there have been some attempts [9-11] to describe the burst atoms using the mean- 
field GP equation, now there seems to be a consensus that they cannot be described 
adequately and satisfactorily using a mean-field approach [12-17]. The GP equation 
is supposed to deal with the zero- or very low-energy condensed phase of atoms and 
has been successfully used to predict the time to collapse, evolution of the collapsing 
condensate as well as its oscillation [7-9,11,12]. However, the GP equation has not been 
fully tested to study the very low-energy (~ nano Kelvin) jet formation [3] when the 
collapse is suddenly stopped before completion by jumping the scattering length to zero 
(noninteracting atoms) or positive (repulsive atoms) values. As the jet atoms are very 
low-energy condensed atoms the mean-field GP equation seems to be suitable for their 
study and we present such a systematic description in this paper. 

In the experiment [3] the initial scattering length ctini|;ig^i (> 0) of a repulsive 
(a > 0) or noninteracting (a = 0) condensate is suddenly jumped to the negative value 
Collapse ^) start the collapse. The condensate then begins to collapse and 
explode. The collapse is then suddenly terminated after an interval of time ^evolve 
by jumping the scattering length from a^^Q^^g^pg^ to Oq^ench (— "^^^ atoms are 
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slowly formed in the radial direction when the collapse is stopped in this fashion. In 
the experiment usually ctq^gj^^^j^ = 0. Sometimes the scattering length is jumped past 
^quench '^expand ~ 250ao to have an expanded condensate to facilitate observation. 
It is emphasized that unlike the emitted uncondensed "hotter" missing and burst atoms 
reported in the experiment [3] the jet atoms form a part of the surviving "colder" 
condensate and hence should be describable by the mean-field GP equation. Yet most 
of the theoretical treatments on the topic [7,8, 11-15] are completely silent about jet 
formation. Saito et al. [9] and Bao et al. [17] present a mean-field description of jet 
formation and Calzetta et al. [16] goes beyond the mean-field model in including the 
effects of quantum field corrections in their description of jet formation. Bao et al. [17] 
employ a fully asymmetric mean-field model and are capable of studying the breakdown 
of axial symmetry in jet formation for experiments of collapse performed in an axially 
symmetric trap [3]. 

Although, the breakdown of axial symmetry and the possible necessity of quantum 
field corrections in the jet formation are interesting topics to be studied carefully, we 
investigate the possibility of explaining the jet formation within an axially- symmetric 
mean-field model. Using the GP equation we study satisfactorily the jet formation for 
the experimental values of the scattering lengths and times which shows that a mean- 
field model describes the essential features of jet formation. The number of jet atoms is in 
agreement with experiment. Further, we extend our study to other values of scattering 
lengths and times and predict the possibility of the formation of jet atoms. Future 
experiments may test these predictions and thus provide a more stringent test for the 
mean-field GP equation. To account for the loss of atoms from the strongly attractive 
collapsing condensate we include an absorptive nonlinear three-body recombination term 
in the GP equation [6]. 

In section 2 we present our mean-field model. In section 3 we present our results 
that we compare with the experiment and other numerical studies. In section 4 we 
present a physical discussion of our findings and some concluding remarks are given in 
section 5. 



2. Nonlinecir Gross-Pitaevskii Equation 



The time-dependent Rose- Einstein condensate wave function \E'(r;r) at position r and 
time T allowing for atomic loss may be described by the following mean-field nonlinear 
GP equation [1] 



X (X2A^|*(r;T)|2 + X3A^2|*(r;T)|^)l*(r;T) = 0. 



(2.1) 



Here m is the mass and the number of atoms in the condensate, g = Anfi^a/m the 
strength of inter-atomic interaction, with a the atomic scattering length. The terms K2 
and K3 denote two-body dipolar and three-body recombination loss-rate coefficients. 
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respectively and include the Bose statistical factors 1/2! and 1/3! needed to describe 
the condensate. The trap potential with cylindrical symmetry may be written as 
V{r) = ^muj^{r'^ + A^z^) where uj is the angular frequency in the radial direction r 
and Xu that in the axial direction z. The normalization condition of the wave function 
is J (ir|\E'(r; r)p = 1. Here we simulate the atom loss via the most important quintic 
three-body term [6-9]. The contribution of the cubic two-body loss term [20] is 
expected to be negligible [6, 9] compared to the three-body term in the present problem 
of the collapsed condensate with large density and will not be considered here. 

In the absence of angular momentum the wave function has the form \&(r;r) = 
il>{r,z;T). Now transforming to dimensionless variables defined by x — y/2r/l, y — 
V2z/l, t — TUJ, I = ^/h/Jrnw), and 



<f){x,y;t) 



^{x,y;t) 



X 



/3 
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(2.3) 



where n — Na/l and ^ = 4X3/ (a^Z^o;). The normahzation condition of the wave function 
becomes 



norm 



= 271 



dx 



dy\ip{x,y;t)\ x 



2^-1 



(2.4) 



For ^3 = 0, Anorm = 1) however, in the presence of loss > 0, A/norm < 1- The 
number of remaining atoms N in the condensate is given by = A^oA/norm, where A^o 
is the initial number of atoms. 

In this study the term will be used for a description of atom loss in the case 
of attractive interaction. Near a Feshbach resonance the variation of K3 vs. scattering 
length a is very rapid and comphcated [18]. From theoretical [19] and experimental [20] 
studies it has been found that for negative a, K3 increases rapidly as |a|", where the 
theoretical study [19] favors n = 2 for smaller values of |a|. In this work we represent 
this variation via a quadratic dependence: K3 ~ a^. This makes the only "parameter" 
^ of the present model a constant for an experimental set up with fixed / and uj and in 
the present study we use a constant ^. In our previous and the present investigation we 
choose ^ or to provide a correct evolution of the number of atoms in the condensate 
during collapse and explosion. The mean-field GP equation is best-suited to make this 
prediction. After a small experimentation it is found that ^ = 2 fits the time evolution 
of the condensate in the experiment of Donley et al. [3] satisfactorily for a wide range of 
variation of initial number of atoms and scattering lengths [7] . This value of ^ is used 
in all simulations reported in this paper. A similar philosophy is Tiscd in choosing the 
value of in [12,17], where the authors reproduced the rate of variation of the number 
of atoms in the collapsing condensate. However, interest in other mean-field [9, 11] and 
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beyond-mean-field [13-15] studies of tlie collapsing condensate was not in reproducing 
the variation of the number of atoms in the collapsing condensate. Consequently, if a 
three-body recombination term is used in these studies other criteria are used in fixing 
the value of K^. 

It is useful to compare this value of 2) with the experimental [20] estimate 
of three-body loss rate of ®^Rb as well as with other values used in the study of 
the experiment by Donley et al. For this we recall that the present value ^ = 2 
with K3 = ^aH^u/A leads to [7, 8] iTg ~ 8 x 10"^^ cmVs at a = -340ao and 

~ 6 X 10^^'' cm^/s at a = — SOao- Santos et al. [11] employed the experimental 
value [20] ~ 7 x 10"^^ cm^/s at a = — 340ao which is very close to the present choice. 
(Santos et al. quote L3 = QK-^. We overlooked this fact in [7,8] and hence misquoted 
there the values of [11,20].) Bao et al. [17] employed ~ 6.75 x 10"^'' cm^/s at 
a = — 30ao in close agreement with the present choice, which leads to a very similar 
time evolution of the number of atoms in the collapsing condensate. Savage et al. [12] 
employed a slightly larger value i^'s ~ 19 x 10~^^ cm^/s at a = — 30ao and also produced 
reasonably similar results for the time evolution of the number of atoms. However, Saito 
et al. [9] and Duine et al. [14] employed a much smaller value (smaller by more than 
an order of magnitude) ~ 2 x 10^^^ cm^/s at a = — 30ao. With this value of K^, 
unlike in the other studies [7, 11, 12, 17], it is not possible to fit the the time evolution 
of the number of atoms in the collapsing condensate. Of these theoretical studies, the 

values used by Santos et al. [11], Savage et al. [12], Bao et al. [17] and the present 
author [7] are consistent with each other and describes well the decay of the collapsing 
condensate. 

3. Numerical Result 

We solve the GP equation ()2.3|) numerically using a time-iteration method based on 
the Crank-Nicholson discretization scheme elaborated in [21]. We discretize the GP 
equation using time step A = 0.001 and space step 0.05 for both x and y spanning x 
from to 15 and y from —40 to 40. This domain of space was sufficient to encompass 
the whole condensate wave function in this study. 

The calculation is performed with the actual parameters of the experiment by 
Donley et al. [3], e. g,, the initial number of atoms, scattering lengths, etc. The 
numerical simulation using ()2.3|) with a nonzero ^(=2) immediately yields the remaining 
number of atoms in the condensate after the jump in scattering length. The remaining 
number of atoms vs. time for ajj^i^ig^i = 7ao, collapse ~ — 30ao, C, = 2, and Nq = 16000 
is in satisfactory agreement with experiment [7]. 

Now we consider the jet formation as in the experiment for these sets of the 
parameters after different evolution times ^gyolve collapsing condensate when 

the scattering length is suddenly changed from ^collapse ~ — 30ao to ctquench ~ ^ 
or to 250ao. First we consider Oq^gj^ch ~ ^- figures 1 (a) and (b) we plot the 
contour plot of the condensate for ^gyolve ~ ^ ™^ ^ respectively, in this case 
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Figure 1. A view of tlie evolution of radial jet at times t = 0, 2ms, 4ms and 5.2 ms 
on a mat of size 16 /im x 16 /im from a contour plot of \(f>{x, y, t)\'^ for ajj^j^jai = 7ao, 
"collapse = -30ao, C = 2, A'^o = 16000, Oquench = 0, and (a) igvolve = 4 ms and (b) 
^evolve = ^ 



at different times t = 0, 2 ms, 4 ms, and 5.2 ms after jumping the scattering lengtfi to 
'^quench ~ ^- ^ prominent radial jet is formed at time t = 5.2 ms after stopping the 
collapse for ^gyolve ~ ^ '^^^ formed slowly after stopping the collapse and is 
more prominent 4 — 6 ms after stopping the collapse. The jet for tg^Qj^g = 8 ms is not 
so spectacular. We also studied the jet formation for ^gyolve ~ ^ ^ 
The jet is much less pronounced for ^gyolve = 2 ms and 10 ms compared to the jet in 
figures 1. 

Next we study the effect of the variation of c^collapse ~ ~30ao on the jet formation. 
For this purpose for the same set of parameters of figures 1 we consider ctgoUapse ~ 
— 6.7ao and o,QQ\\g^p^Q = — 250ao corresponding to smaller and larger attraction, in figures 
2 (a) and (b), respectively. First, we consider the case ^collapse ~ — 6.7ao in figure 2 
(a). In this case the final attraction is weaker and the collapse is less dramatic. The 
collapse and the decay of atoms do not start until ^gyolve ~ ^ -^^^ ^evolve ^ ^ ^ 
broad jet is formed after about 12 ms of disruption of collapse. After 4 ms of disruption 
of collapse there is almost no visible jet. The jets are broader and less prominent for 
larger values of ^gyolve' Next we consider the highly attractive case OgQjjg^pgg = — 250ao 
in figure 2 (b). In this case due to a very large attraction the collapse and the decay of 
atoms start at a small value of ^gyolve '^^'^^^ zero. Hence a reasonable jet is formed for 
^evolve ~ ^ small times after stopping the collapse. Because of large attraction 

in this case the collapse is over for a smaller value of ^gyolve' Hence there is almost no 
jet formation for ^gyolve — ^ However, the nature of jets in each case of figures 1 
and 2 is distinct. 
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Figure 2. A view of tlie evolution of radial jet on a mat of size 16 /im x 16 fim from 
a contour plot of |(/>(a;, y, f)p for ajj^j^jaj = Toq, ^ = 2, A'q = 16000, ^quench ~ ^'^'^ 
(a) ^collapse ~ — 6.7ao, ^evolve ~ ^ times i = 0, 4ms, 8ms and 12 ms and (b) 

•^collapse ~ — 250ao, ^evolve ~ ^ times t = 0, 2ms, 4ms and 5.2 ms. 



Donley et al. [3] also considered expanding the condensate before observing the 
jet by jumping the scattering length to a large positive value ^expand ~ 250ao after 
disrupting the collapse. This procedure expands the size of the condensate so that it 
might be easier (or the only way) to observe the jets in the laboratory. However, we 
find that in all cases the jet is much less pronounced after this expansion. This is 
illustrated in figures 3 (a) and (b) where we plot the jet formation corresponding to the 
cases reported in figures 1 (a) and (b), respectively, after expanding the condensate to 
'^expand ~ 250ao as in the experimental result reported in figure 5 of [3]. In plots of 
figures 3 (a) and (b) the condensate is of larger size than in the corresponding plots of 
figures 1 (a) and (b). In figure 3 (a) the jet is almost destroyed. In figure 3 (b) the jet 
appears but it is wider due to expansion. We also expanded to Og^pg^j^^ = 250ao the 
jets for ^collapse ~ — 6.7ao for different ^evolve' J^^^ were almost destroyed in those 
cases. After expanding to ^expand ~ 250ao the jets for o-^oHapse ~ — 250ao became 
wider due to expansion but remained visible. 

In addition, we studied jet formation for different values of ajni^^iai in place of 
'^initial ~ "^^^ ^"^^ ihaX the scenario remains very similar independent of the initial 
scattering length. However, the number of particles in the jet gives a quantitative 
measure of jet formation and in the following we make a study of the number of atoms 
in the jet in different cases. 

Next we calculated the number of jet atoms in different cases by integrating the 
wave function over the relevant region where the jet is formed. The normalization 
condition ()2.4p gives the total number of atoms in the condensate via Ai'oA/norm- After 
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Figure 3. A view of the evolution of radial jet at times t = 0, 2ms, 4ms and 5.2 ms 
on a mat of size 16 /im x 16 /im from a contour plot of \(f>{x, y, t)\'^ for ajj^j^ja,! = '''^Oi 
"collapse = -30ao, f = 2, iVo = 16000, a^-j^pand = 250ao, and (a) igyoiye = 4 ms and 
(b) ^evolve = ^ ms. 



an examination of figures 1 and 2 or other relevant jet figures we separate the condensate 
at 5.2 ms (as in the experiment) according to x values into the central part and jet. 
The X integral in ()2.4|) is then separated into the central part and jet in each case. 
The jet part (outer x values) of the integral p.4|l multiplied by A^^o gives the number 
of jet atoms. A similar procedure has been used in [17] to calculate the number of jet 
atoms. However, in actual experiment to see jet or any other phenomenon, the harmonic 
trap has to be removed and the condensate allowed to expand and photographed. This 
enlarges the condensate to be photographed without presumably losing its actual shape 
and characteristics. Vortices and dark and bright solitons photographed in this fashion 
give the true picture of the condensate before free expansion. In numerical simulation, 
on the other hand, it is possible and much easier to count the jet atoms more accurately 
without any expansion. Assuming that there is not much experimental error in counting 
the jet atoms after free expansion we attempt to compare the two in the following. 

For a fixed A^o = 16000 and ^initial — '^^o we calculate the number of atoms in the 
jet for different evolution time ^gyolve '^collapse" none of the cases an expansion 
to Ogjjpg^j^^ was applied. The results are plotted in figures 4. In figure 4 (a) we plot the 
variation of the number in jet vs. ^evolve '^collapse ~ — 30ao and in figure 4 (b) we 
plot the variation of the number in jet vs. Ic^collapsel/'^o '^evolve ~ ^ figure 
4 (a) we also plot the experimental result for the number of atoms in jet for the same 
values of the parameters and find the agreement of our calculation with experiment to 
be quite satisfactory. It is noted that an expansion to ^expand applied in the 

data reported in figure 6 (a) of [3,22]. 
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Figure 4. Number of atoms in jet 5.2 ms after jumping tlie scattering length 
*° "^quench = ^ = 16000, a-^i^[t'm\ = "^"o, £, = 2 vs. (a) igvolve 

'^collapse ^ — 30ao and (b) I acoUapse I /"^o ^evolve ~ ^ ™' ^^^^ ~ experiment 
of Donley et al. [3], black dashed line — mean field model of Bao et al. [17], blue dotted 
line — average over theoretical results of Calzetta et al. [16], red full line — present 
result. 

In addition, in figure 4 (a) we compare the present results with the theoretical 
calculations by Calzetta et al. [16] and Bao et al. [17]. The calculation by Bao et al. [17] is 
essentially based on mean-field GP equations as in this study and they correctly identify 
the jet atoms as being a part of the condensate. They employ a fully asymmetric mean- 
field model in their description of jet formation. Yet the present result for the number 
of atoms in the jet is larger than theirs and in better agreement with experiment [3]. 
The reason for this discrepancy is unknown. 

The calculation by Calzetta et al. [16] uses a theoretical model beyond mean-field 
taking into consideration quantum fiuctuations. The number of jet atoms of Calzetta et 
al. [16] is in good agreement with the experiment [3] as well as with the present study. 
In figure 4 (a) we plot the average over the number of jet atoms obtained by Calzetta 
et al. As the physical inputs and the dynamics of the study of Calzetta et al. and 
the present study are quite distinct, it is difficult to compare the two and conclude 
about the effect of quantum fluctuations on jet formation. The effect of quantum 
fluctuations could turn out to be significant in various aspects of the experiment of 
the collapse of a BEC including the formation of jet atoms. For example, they are 
of utmost relevance in the molecule formation in a collapsing BEC near a Feshbach 
resonance and in subsequent rapid atom- molecule oscillation as observed recently [23]. 
Also, the burst atoms cannot be properly described by a mean-field model and quantum 
corrections could be significant [16,17]. Further studies are needed to identify clearly 
the effect of quantum fiuctuations on jet formation. 

The investigation by Saito et al. [9] is very similar to this study in applying an 
axially-symmetric mean-field model and producing the essentials of jet formation. They 
also provide a physical explanation of jet formation. In the collapsing condensate two 
distinct spikes are formed in the condensate wave function along the axial direction as 
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the atomic interaction is changed from repulsive to attractive. These two spikes act as 
sources of matter waves and jet is the interference pattern of matter waves from these 
two sources [9]. However, we feel that small local spike(s) in the wave function expand 
in the radial direction in the form of a jet when the collapse is stopped by removing the 
atomic attraction and the force on the spike (s) is suddenly changed [3]. 



4. Discussion 



In this section we give a physical explanation of the different types of jet formation 
noted in last section. The jet formation is more dramatic when the vigorous collapse of 
a condensate is suddenly stopped by turning the attractive condensate noninteracting 
(^quench ~ ^) repulsive (c^cxpand ~ 250ao). In the strongly collapsing condensate 
local radial spikes are formed during particle loss as can be seen from a plot of the 
numerically calculated wave function [7] and in experiment [3] . During particle loss the 
top of the spikes are torn and ejected out and new spikes are formed until the explosion 
and particle loss are over. There is a balance between central atomic attractive force and 
the outward kinetic pressure. If the attractive force is now suddenly removed by stopping 
the collapse by applying Oquench ~ ^' highly collapsed condensate expands due to 
kinetic pressure, becomes larger and the recombination of atoms is greatly reduced. 
Consequently, the spikes expand and develop into a prominent jet [3] for Oq^gj^^^j^ = 
as in figure 1 (a). 

However, if the condensate is expanded further by applying Ogxpand ~ 250ao, in 
the cases studied, the spike as well as the condensate expand so much that the prominent 
jet becomes in general more diffuse as in figure 3 (b) or completely destroyed as in figure 

3 (a). If the attractive condensate is allowed to collapse for sufficiently long time, the 
explosion stops eventually and a relatively cold remnant condensate is formed. At that 
stage there would be almost no prominent spikes in the wave function and no jet could 
be formed by applying Oq^ench ~ ^- With the increase of evolution time of a collapsing 
condensate the jet becomes less prominent as can be seen in figures 1 (a), (b) and 4 (a) 
and eventually disappears for large ^gyolvc However, the collapse and decay of particles 
start after a finite ^pyolve before which the jet formation is practically absent as in figure 

4 (a). If the condensate is weakly attractive as in figure 2 (a) the collapse is also weak 
and the spikes are less pronounced. Consequently, upon stopping the collapse a wide jet 
is formed after a longer interval of time as in figure 2 (a). As I "^collapse I increased at 
a fixed ^evolve = 4 ms as in figure 4 (b) one gradually passes from a strongly collapsing 
condensate to a relatively cold remnant, that is from a region of prominent jet formation 
to a region with less prominent jet. For smaller Ic^collapsel collapse is weaker and 
jet formation is absent as in figure 4 (b). 

In the actual experiment [3] the jet is found not to possess axial symmetry. In the 
present study we use an axially symmetric model to describe the essential features of the 
jet. Hence, although an axially symmetric model is enough for a qualitative description 
of the jet, a full three-dimensional model might be necessary for its complete quantitative 
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description. 

5. Conclusion 

In conclusion, we have employed a numerical simulation based on the accurate solution 
[21] of the mean- field Gross-Pitaevskii equation with a cylindrical trap to study the jet 
formation as observed in the recent experiment of an attractive collapsing condensate 
by Donley et al. [3]. In the GP equation we include a quintic three-body nonlinear 
recombination loss term [6] that accounts for the decay of the strongly attractive 
condensate. The result of the present simulation is in good agreement with the 
experimental result for jet formation [3]. We also compare the present result with two 
other recent theoretical calculations of jet formation [16,17]. Of the different aspects 
of the experiment by Donley et al. the dynamics of relatively hot emitted burst and 
missing atoms seems to be beyond mean-field treatment [16,17]. However, the various 
properties of the cold residual condensate including jet formation seem to be describable 
by mean-field models. In fact, many features of the experiment by by Donley et al. [3], 
specially the detailed behavior of the surviving remnant condensate [7, 9, 17] and jet 
formation, have been understood by introducing the rather conventional three-body 
recombination loss in the standard mean-field GP equation, with a loss rate compatible 
with other studies [7, 12, 17, 18, 20]. 
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